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Abstract 

A general theory is presented for the treatment of decoherence of a multilevel quantum system (with 
many degrees of freedom) interacting with multi-bath reservoir and driven by ac fields. In this approach, 
the system is described by a reduced density operator and the multi-bath reservoir is characterized by a 
number of spectral densities. The reduced density operator is governed by the master equation in which the 
effect of ac driving fields and the leakage to non-computational states are included. The theory is applied 
to the study of decoherence of a two-dimensional (2D) SQUID flux qubit coupled to control and readout 
circuits. The predicted results are in very good agreement with available experimental results in the absence 
of driving fields and with the analytic results of a dissipative two-level system in the presence of weak 
driving fields. The relaxation and decoherence times versus the parameters and temperature of the control 
and readout circuits are also explored in details to facilitate the optimization of the 2D SQUID qubit. 
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I. INTRODUCTION 



In the past few years, solid-state qubits based on superconducting devices are of particular 
interest for quantum computation because of their advantages of large-scale integration, flexibil- 
ity in design, and easy connection to conventional electronic circuits^. Significant progress has 
been made on physical implementation of quantum computation based on superconducting qubits. 
Quantum coherence has been successfully demonstrated in a variety of superconducting single- 
qubit systema^'Sii^i^'^i^i^i^ and coupled two-qubit systems.— >iiii^ii^ However, all the superconduct- 
ing qubits demonstrated in the experiments so far have relatively short coherence time and high 
probability of gate errors ^i^iii^t^i^ i ^ One of the causes of these problems is the intrinsic gate 
error resulting from leakage to non-computational states due to the typical multilevel structures 
of superconducting qubits.— This kind of gate error has been explore d^^'^^'^^ and can be mini- 
mized by using appropriate working parameters for the qubits with their device parameters given 
in priori^ Another cause is the extrinsic gate error arising from coupling between the qubits and 
environment resulting in decoherence such as relaxation and dephasing.— Due to the unavoidable 
coupling with environment, the superconducting qubits always suffer from such kind of extrinsic 
gate error. Thus the environment-induced decoherence is one of the main obstacles to the practical 
application of superconducting qubits in quantum computation.— i^ i^^i^^i^*^ 

The typical environment in superconducting qubits is electronic circuits used for con- 
trol and readout of the qubits. Although the decoherence of superconducting qubits 
induced by such kind of environment has been extensively investigated both theoreti- 
callyldld§^^^M.^^^^^^^ML2MMi^ and expeTimentally^^^^^^^M^^^^^^, almost 
all the investigations so far are for the qubits in the absence of driving fields (free decay). However, 
in superconducting-qubit based quantum computation, ac fields (e.g., microwave fields) are usu- 
ally used to manipulate the qubit's states^'^'^'^ ^^^^^^^^^^'^^^^^'^^ Due to the coupling of the qubits 
with the fields the effect of driving fields and leakage to non-computational states may be quite 
large depending on the field strengtb^^'i^'i^ and thus the dynamics of driven qubits may be quite 
different from that of qubits in free decay. Recent experiment's shows that the decoherence time 
of a superconducting qubit is significantly increased in the presence of a resonantly ac driving field. 
Thus a through investigation of decoherence of realistic superconducting qubits needs to include 
effect of driving fields and leakage to non-computational states.— 

Furthermore, building a practical quantum computer requires simultaneous operations of a 
large number of multiqubit gates in a coupled multiqubit system.— i^i^ On one hand, the coupled 
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multiqubit system may have many degrees of freedom, on the other hands, due to complexity 
and diversity of the superconducting circuits, the environment may be a multibath one. In this 
case, the coupled multiqubit system may interact with the environment through all the degrees 
of freedom simultaneously, resulting in significantly different decoherence from that interacting 
with the environment through one degree of freedom of single qubits. Therefore, to investigate the 
decoherence of coupled superconducting qubits in realistic gate operations of quantum computation, 
one needs an approach for a quantum system of many degrees of freedom interacting with a 
multibath reservoir and driven by ac fields. 

In this paper, we propose a general theory for the treatment of decoheerence of a multilevel 
quantum system of many degrees of freedom interacting with a multibath reservoir and driven 
by ac fields. In this theory, the multibath reservoir are characterized by a number of spectral 
densities.— 1^1^ For superconducting qubits coupled to electronic circuits, the spectral densities at 
finite temperature can be calculated in terms of recently proposed approaches^ii^ together with 
the quantum fluctuation-dissipation theorem.— The system is described by a reduced density 
operator. The reduced density operator is governed by the master equation in which the effect of 
driving fields and leakage due to the driving field and reservoir are included. This theory is used 
to simulate the dynamic process of a two-dimensional (2D) superconducting quantum interface 
device (SQUID) flux qubit coupled to control and readout circuits in the absence of driving fields 
(free decay). The results are in very good agreement with the available experimental results. It is 
also applied to investigate the effect of driving fields on the decoherence of the 2D SQUID qubit 
coupled to the control and readout circuits and driven by a resonant microwave field. The results 
agree well with the analytical results of the dissipative two-level system in the case of weak driving 
fields. To optimize the 2D SQUID qubit, the changes of relaxation and decoherence times versus 
the parameters and temperature of the control and readout circuits are explored. 

II. MASTER EQUATION OF A DRIVEN OPEN QUANTUM SYSTEM 
A. General form of master equation for a driven open quantum system 

In general, an open quantum system is described by a generalized master equation of non- 
Markovian process.— i^i^i^ However, in the case of weak damping and weak driving fields, which 
is the case for most of superconducting (charge, flux, and phase) qubits, the generalized non- 
Markovian master equation is equivalent to the Markovian master equation.— Thus we will 
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present the dissipative theory for Markovian process. 

Let us consider a global system consisting of a quantum system surrounded by a reservoir 
and driven by an ac field. If the interaction of the driving field and reservoir is neglected, the 
Hamiltonian of the global system can be written as 

Hg = Hs + Hr + Hi + Hp, (1) 

where, Hs is the Hamiltonian of the quantum system, Hji is the Hamiltonian of the reservoir, Hi 
is the interaction of the quantum system and reservoir, and Hp is the interaction of the quantum 
system and driving field. Obviously Hr commutes with both Hs and Hp. For the case of weak 
damping and weak driving field considered here, the interactions Hj and Hp are proportional to 
the system's coordinate operators linearly and thus they also commute with each other.— 

In Schrodinger picture, the motion of the global system is described by Liouville-von Neumann 
equation^ 

^ = [HG{tU{t)] = -tC{t)r,{tl (2) 

where, r]{t) is the density operator of the global system and C{t) is the Liouville superoperator 
corresponding to the Hamiltonian of the global system defined by 

C =i [Hg, ]= Cs + Cr + Ci + Cp. (3) 

Here, Cq = [Hq, ] /h is the Liouville superoperator corresponding to the Hamiltonian Hq for q = S, 
R, I, and F. The density operator of the global system satisfies Trs,Rf]{t) = 1, where Tis,r is the 
trace over both the quantum system and reservoir. 

In the global system, the quantum system performs a dissipative process due to the coupling with 
reservoir. This process can be characterized by a reduced density operator p{t). It is calculated 
by tracing r]{t) over the reservoir 

p{t) = TvRijit). (4) 

In the case of weak field and weak damping, p{t) is governed by the master equation which, in 
Schrodinger picture, is given by (see Appendix A for details) 

^ = -z [Cs + Cp{t)] p{t) + Vj{t)p{t), (5) 

where, P/ {t) is the dissipation superoperator given by Eq. ()A13p . On the right-hand side of Eq. 
([5]), the first term describes the coherent dynamics of the unperturbed quantum system, the second 



term represents the pumping dynamics of the driving field, and the third term which is called the 
dissipator— characterizes the dissipative dynamics due to the interaction of the quantum system 
and reservoir. The dissipator includes all the information of dissipative dynamics;^ Note that 
Eq. ([5]) is determined only by the operators at present and no longer depends on the operators 
in the past. Thus it describes a Markovian process of the quantum system^-, ^he Markovian 
process is irreversible because the energy transferred to the reservoir can no longer return to the 
quantum system completely. In the case of weak damping and weak driving fields the Markovian 
approximation is equivalent to the non-Markovian approach.— 

B. Master equation of a driven quantum system of many degrees of freedom interacting 
with a multibath reservoir 

Let us consider a quantum system of many degrees of freedom encompassed by a multibath 
reservoir and driven by an ac field. The multibath can be separated into a number of independent 
baths. They have different properties and interact with the system via different degrees of freedom. 
If interactions between the baths are incorporated into the Hamiltonian of the reservoir, the total 
interaction Hj of the system and the multibath reservoir is given 



where, is the operator of the ^uth bath of the reservoir, X^{t) is the system operator cor- 

responding to the /ith bath, and the sum to /i is over all the baths of the reservoir and all the 
degrees of freedom of the system. In general, the interaction Hi{t) is a Hermitian operator but 
Xfj_ and may be non-Hermitian. In the case of weak damping, is reduced to the system 

coordinate operator and <?^(t) is reduced to the fluctuation force of the /xth bath.— For a driven 
quantum system of many degrees of freedom interacting with a multibath reservoir, the reduced 
density operator of the system still satisfies the maser equation given by Eq. ([5]). Using Eqs. ([3]), 
dl]), dM]), (PT3l) . and ([XTi]) . the dissipator Vi{t)p{t) in Eq. jS]) is now given by 




(6) 




(7) 
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where, is the representation of in the interaction picture given by 



(-r) = exp ( '-Hst\ X^ {t + r) exp [-^Hst ] , (8) 



and J'fj_i, (r) and {—t) are the two kinds of correlation functions of the fluctuation forces between 
the //th and z^th baths. The correlation function J^i, (r) is defined by 

J'^^(t) = Tr^j[?^(t)a(i?)?^(t-r)] 

= (?/.(t)?.(t-r)), (9) 

where, ?^ is the representation of in the interaction picture defined by 

?M (r) = exp ^^Hrt^ (t) exp ^--^i^/jr^ . (10) 

The correlation function J'^^ {—t) is the complex transpose of ^^j, (r). In general, since and 
— t) may not be Hermitian operators and/or they do not commute with each other, jT'^j^ (r) is 
a non-Hermitian matrix and J'u^ {—'t) 7^ -Jfiu (t). 

If the //th bath does not interact with the z^th bath, the correlation between the two baths is 
zero. In this case, the correlation function J^i, (r) can be written as 

Jmv (r) = J^f, (r) (11) 

where, (r) is the autocorrelation function. The autocorrelation function (r) represents the 
correlation of the bath state at time t — r to the bath state at time t, while the autocorrelation 
function J'/l^ {—'t) represents the correlation of the bath state at time t to the bath state at time 
t — T. In general, J^/l^ {—'t) 7^ •Jfifi {t), which signifies the irreversibility in time of the correlation 
function. From Eq. ([9]), the irreversibility is attributed to the non-Hermitian property and/or 
noncommutable property of the reservoir operators. 

From Eq. ([7]), the effect of reservoir on the system is represented by the correlation function 
J^{t). It is difficult to calculate the correlation function directly from the fluctuation forces of the 
reservoir. In reality, the effect of reservoir on the system can be equivalently characterized by a 
real and measurable spectral density of the reservoir— i^i^. The relation between the correlation 
function J^{t) and the spectral density [lo) is given by the Wiener-Khintchine theorem^ 

J{t) = — / (cj)exp(fcjt)(ia;, (12) 

and 

/ + 00 
J{t) exp {-iujt) dt. (13) 
-oo 



For solid-state qubits, particularly for the superconducting qubits consisting of electronic cir- 
cuits, the spectral density at zero temperature can be calculated in terms of the recently pro- 
posed approaches j^i^ The spectral density at finite temperature can be obtained from the zero- 
temperature spectral density by means of the quantum fluctuation-dissipation theorem^i^. 



III. REPRESENTATION OF MASTER EQUATION IN HILBERT SPACE 
A. Representation in an arbitrary Hilbert space 

Suppose {\n)} is a complete basis set of the system. In the Hilbert space spanned by {|n)}, the 
reduced density operator p is represented by a density matrix with matrix elements given by 

Pmnit) = {m\pit)\n) ■ (14) 

The diagonal matrix element /o^m off-diagonal matrix element (m 7^ n) are the population 
of state \m) and coherence of states \m) and respectively. From Eqs. ([5]) and ([7]), the density 
matrix element (t) satisfies the master equation 



- / J i~ rnn,rri- rf ■ — mn,m'n' \ rm'n' 
m'n' 

+ '^iin,m'n'Pm'n'^ (15) 
m'n' 

where, Cf^n,m'n' ^mn,m'n' matrix elements of Cs and Cp, respectively, and T^mn,m'n' 

is the dissipation matrix element of the dissipation superoperator Vf describing the effect of the 
reservoir on the system. They are given by 

^mn,m'n' ~ [^mm'^n'n ~ H~,^Smm'~\ ; (16) 

and 

ft 



'^mn,m'n' — ^ / \_^mm' i'^) ^n'n i''') '^fJ-'^ i''') 

~Snn' J2 ^-fc ^km' (^) (-^) 
k 

-6mm'Y.^n'k M Xj^^ (t) J^, (t) 
k 

+ (r) (t) J^V (-^)] dr, (17) 

where, = {m\HQ\n) for = S and F, and (t) = {m\Xa {t)\n) and (t) = 



m 



n) for a = ji and v. 



B. Representation in the Hilbert space of Hamiltonian eigenstates 



Suppose En and |n) are respectively the eigenvalue and eigenfunction of Hs obtained by solving 
the eigenvalue equation Hs \n) = En\n). In the Hilbert space spanned by the eigenstates {|n)}, 
the density matrix element and the master equation are still given by Eq. (jl4p and Eq. (jlSp . 
respectively. But now the matrix elements and are simplified to 

and 

where, uJmn = {Em — En) /h and H^nit) = \HF{t)\ n). The matrix element T^mnm'n' 
p!7|) is also simplified to 

-^nw Yl Kk it) Xlm' it') e'^-'^'Jl i-r) 

k 

-5mm' K'k it') Kn it) e^"^-' (^) 
k 

+ X^n^.^' K'n it) e^"™'™ ("T)] dr , (20) 

where, t' = t — t. 

If the Hamiltonian eigenfunctions {|n)} are real, one has -ff^n ~ ^nm Hermitian operator 

Hp{t). From Eqs. ([18|) and (fT9]) . the matrix elements i2®^ for B = S and F have the following 
symmetric relations 

'~'nm,n'm' '~'mn,m'n'j \ J 

and 

'~"m'n' ,mn '~'mn,m'n'' \ ) 

From Eq. ([2T]) one has ll%^m,nn = 0- 

If is time-independent, the matrix elements Dlnn m'n' ^^^^ further simplified. Introduce a 
function F^^{uj,t) by 

F^,(u;,t)= f e''^^J^,{T)dT. (23) 



Applying Eq. to Eq. ([20]) one obtains 

k 

-^mm' ^ Xn'k^kn^f^''i^kn' ,t) 

k 

+ X^m'KnFU^m^'n'^t)] , (24) 

where F^,^{u!,t) is the conjugate transpose of F^iy{uj,t). If the interaction between any two baths 
is neghgible, T^mnm'n' '^^^ further simphfied by using Eqs. (fTTj) and ([^5|) . In this case, the 
correlation terms between different baths in Eq. (j24p vanishes and thus the dissipation matrix 
element of the multibath reservoir is a sum of individual dissipation matrix element of each bath 
of the reservoir. 

IF time-dependent part of is separable from system operators, can be expressed in a 
general form of Fourier series as 

^MW = E^Ae^PK*)' (25) 
A 

where, X\ is a time-independent operator that only depends on system operators. 
Substituting Eq. ([25]) into Eq. ([20]) and using Eq. ([231), one has 



'^mn,m'n' ~ »-2 ^-^P 



fJ-,!^ A1A2 



M t^Ai + ^Aa I t 



^ [^Ai, mm' XA2,n'n-^Mi' (^A2,nn"* 



-^mm' X^XA2,n'fcXAi,fcn-^M'^ (^A2,A:n"*) 
k 

+ XA2,mm'XAi,n'n^Mi^ i'^A^mm"*) ' (^6) 



where, <^'^xfnn' ~ ^ '^A2 Xxmn ~ I^aI matrix element of the operator x"- 

Furthermore, if Xx is a time-independent operator that does not depend on system operators 
we have XA,mn = ("'-IXaI"') = XA'^mn- In this case, from Eq. ([26]) one has T^im.m'n' = 0- Thus 
if the interaction of the system and reservoir does not depend on system operators, the reservoir 
does not have any effect on the system and thus they are decoupled. 
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IV. DISSIPATION OF A DRIVEN QUANTUM SYSTEM DUE TO A THERMAL BATH 



One of the most popular and important reservoirs is a thermal bath. Due to the interaction 
with the thermal bath, the quantum system transits from one thermodynamic equilibrium state to 
another equilibrium state.^^- During this process, the spontaneous decay and stimulated transition 
follow the detailed balance principle. As has been demonstrated, if the interactions between the 
independent baths are neglected the dissipation matrix element of the multibath reservoir is a sum 
of individual dissipation matrix element of each bath. Thus in this section we will present the 
dissipative theory for a driven quantum system interacting with a thermal bath. 



A. Lamb shift matrix and damping rate matrix 

Suppose the driven quantum system is surrounded only by a sufficiently large thermal bath 
and the system operator X such as the system's canonical coordinate does not depend on time 
explicitly. From Eq. (j24p the dissipation matrix element due to the thermal bath is now given by 

{ k 

+X„im'Xn'n [F{^nn','t) + F* {id rnm' , t)] 

— ^mm' y^^Xn'kXknFiuJkn'it) > , (27) 
k J 

where, the superscripts fi and v are omitted for simplicity, F{uj, t) is given by Eq. ()23p . and F*{uj, t) 
is the complex conjugate of F{uj^t). 

If the spectral density of the thermal bath at temperature T is denoted by (a;), the autocor- 
relation function J (t) can be calculated from (uj) by Eq. (fT2]) . Substituting Eq. (fT2]) into Eq. 
(f23ll we obtain 

F{Lo,t) = FR{oj,t)+iFi{u;,t), (28) 
where, the real part and imaginary part of F(uj,t), Fji{uj,t) and Fj{uJ,t), are given by 

-j^ r+oo 

Fe{uj, t) = — J, [oj') A© [J + uj, t) duj\ (29) 

'^'^ J —oo 

where, Q = R and /, and (lo, t) = sin (ut) /u and A/ (w, t) = [1 — cos (ujt)] juj are the real part 
and imaginary part of A (a;, t) given by 

A (w, I) = /* e^'^Mr. (30) 
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Substituting Eq. ([28]) into Eq. ([27]) we obtain 

where, T^mn,m'n' ■, which we call damping rate matrix element, is the rate matrix element related to 
the change of density matrix from the value of p^'n' to the value of Pmn~ 

+Xmm'Xn'n [FR{uJnn',t) + Fji{uJrnm',t)] 

— Smm' ^ XnikXknFR{iOkn',t) > , (32) 

k } 

and Bmn,m'n' is the Lamb shift matrix element which leads to Lamb shifts of unperturbed energy 
levela^i^ 

{ k 

— ^mm' Xn'kXknFl (cJfcn' , t) > . (33) 

k J 

They satisfy the symmetric relations 



and 



From Eq. ([35]) Bmm,nn (t) = 0. 



B. Steady Lamb shift matrix and damping rate matrix at a long time limit 

In general, the damping rate matrix and Lamb shift matrix are time-dependent because both 
Fji{uj,t) and Fj(uj,t) depend on time. However, in the investigation of decoherence of a qubit, 
one is only interested in the behavior of the qubit after a sufficiently long time. In this case, the 
damping rate matrix and Lamb shift matrix can be well approximated by steady ones. 

In fact, the autocorrelation function of a thermal bath impacts the system only in some time 
interval tc which is called the correlation time^^. As long as the upper limit in the integral of Eq. 
([30|) t » tc it may be extended to the infinity with very little error. 
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Applying the integral 



/ 

Jo 



■oo 



1 

OJ 



(36) 



where V represents the Cauchy principal value of the integral, to Eq. ([50]) one has 



Xr (cj, oo) = 7r5 (tj) , 
Xi {oj, co) = V—. 



(37) 



(38) 



Substituting Eqs. (j37|) and ([38]) into Eq. ([29|) . one obtains 



fuiu) = Fr{uj, oo) = -J^{-uj) , 



(39) 



fiiio) = Fj{u;, oo) = —V 



1 




duj'. 



(40) 



In Eqs. (f32|) and (f33l) . replacing Fr with fa and F/ with //, we obtain steady damping rate 
matrix element Rmn,m'n' and Lamb shift matrix element Bmn,m'n'- They are given by 



respectively. The Rmn,m'n' and Bmn,m'n' have the same symmetric relations as TZmn,m'n' and 
Smn,m'n' given by Eqs. ([MI) and ([35]). 

It is shown from Eq. (j3ip that the thermal bath affects the system via the damping rate matrix 
and Lamb shift matrix. The damping rate matrix elements represent decay rates in special cases. 
For instance, Rmm,mm characterizes the rate of population change in the state m, Rnn,mm describes 
the rate of population transfer from the states m to n, and Rmn.mn represents the dephasing 
(coherence decay between the states m to n) rate of the off-diagonal elements of p^^. In contrast, 
the Lamb shift matrix element such as Bmn,mm represents the well-known Lamb shift of the state 
m induced by the thermal bath. Hence, the damping rate matrix makes the system relax and 




(41) 



and 




(42) 
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decohere while the Lamb shift matrix makes the energy levels shift. In addition, the symmetries 
of the Lamb shift matrix given by Eq. (|35p are the same as those of the matrix for driving fields 
given by Eq. (j2ip . In particular, when substituting Eq. (jSip into Eq. (|15p the Lamb shift matrix 
can be incorporated into the matrix ^i^, and leads to the renormalization of the quantum 
system Hamiltonian.^^ Thus the effect of the Lamb shift matrix is analogous to an extra field. In 
the case of weak damping, the Lamb shift matrix is very small compared to the driving field and 
thus is neglected hereafter. 



C. Spectral density, spontaneous decay, stimulated transition, and detailed balance 

The fluctuation of the thermal bath at temperature T is characterized by a spectral density. 
The spectral density is computed from the real part of the frequency-dependent damping coef- 
ficient 7^ (uj) in terms of the quantum fiuctuation-dissipation theorem (see the Appendix B for 
details).—'^ For the interaction given by Eq. ([6]), the fiuctuation force ^ (t) is proportional to the 
fiuctuation force <; (t) by 



e (t) = A? it) , 



(43) 



where, 



A 



dX 

dq ' 



(44) 



and q is the coordinate operator of the system. For the weak damping considered here the in- 
teraction given by Eq. (I6|) is a linear function of system coordinate operator and thus A is a 
constant. 

For the thermal bath at temperature T, the spectral density (uj) of the fiuctuation force ^ (t) 
is given by Eq. ()B3p . From Eqs. ()B3P and (03]), the spectral density (lu) of the fiuctuation force 
q (t) is given by 



(45) 



A2 A2 

Substituting Eq. (I45p into Eq. (jl2p . we obtain the autocorrelation function of the fiuctuation force 

Mh 



Jit) 



27rA2 



uJjR (uj) 



X exp {iujt) dw. 



(46) 
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In general, the autocorrelation function given by Eq. (|46|) is complex and irreversible since the 
integrand does not have a definite symmetry with respect to to. This result is totally different 
from that for a classical system, for which the autocorrelation function is real and reversible. As 
has been demonstrated, the irreversibility of the autocorrelation function results from the non- 
Hermitian property and noncommutable property of the thermal bath operators. 

Substituting Eqs. (j45p and (j39p into Eq. ()4ip . we obtain for m ^ n the rate of population 
transfer from the state n to m 



1 — coth 



(47) 



\2kBT J 

From Eqs. (j4ip and (j47p we also obtain the rate of population change in the state n 

Rnn,nn — ^ ^ R'mm,nn- (^^) 
m(^n) 

This indicates that Rnn,nn is a sum of the rates of population transfer from the state n to all the 
other states m. This result is identical with that of atomic systems.— 

1. Spontaneous decay 

Let us consider a pair of states denoted by m and n. Suppose the state m has lower energy than 
the state n. The transition frequency is ojnm = {En — Em) /?i > 0. If the system is in the higher 
energy state initially, it can decay without presence of any external field due to the stimulation of 
thermal fluctuation of the thermal bath. This process is the so-called spontaneous decay. During 
this process the system transits from the higher energy state n to the lower energy state m with 
emission of radiation. From Eq. the spontaneous decay rate P;^^ is given by 

^ mn — ^nnilR\^nm) 



1 + coth ' ' 



(49) 



\2kBT J 

It is proportional to the transition matrix element X^n, energy level spacing u>nm, and damping 
coefficient ^^.i^nm)- It also depends on the temperature via the factor coth. (huj nm/ 2k bT). 

2. Stimulated transition 

If the system is in the lower energy state initially it may occur an inverse process with respect 
to the spontaneous decay under the stimulation of thermal fluctuation, making the system transit 
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from the lower energy state m to the higher energy state n. This process is the so-called stimulated 
transition. From Eq. (jTr|) the stimulated transition rate T^J^ is given by 



1 — coth 



WbT) 



(50) 



3. Detailed balance 

In general, (w) is an even function of uj, i.e., ^^^{—lj) = 7^(0;). Using Eqs. (jl9]) and ([50] 
one has 

-pSP /TP p 

mn / -'^n J^m 



Thus the ratio of the spontaneous decay rate to the stimulated transition rate satisfies the detailed 
balance principle. 

V. DECOHERENCE OF A DISSIPATIVE TWO-LEVEL SYSTEM 

An ideal qubit is a two-level system and an ideal qubit interacting with a thermal bath is 
equivalent to a dissipative two- level system (DTLS). For a DTLS, the master equation (jl5p can be 
replaced by the Bloch equation.^'' In the rotating- wave approximation (RWA), analytical solutions 
of the Bloch equation can be obtained, from which analytical expressions of characteristic (relax- 
ation, decoherence, and dephasing) times of the DTLS can be derived. In this section, we present 
some important analytical results without derivation. The details about these analytical results 
can be found in our previous paper.— 

A. Relaxation and decoherence times of the DTLS in the absence of driving fields 

In the absence of driving fields, the relaxation and decoherence times of the DTLS are given by 

Ti = K^\ (52) 

and 

T2 = (53) 
respectively, where, ki is the relaxation rate of the DTLS in free decay given by 

= -^22,11 + -Rll,22, (54) 
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and K2 is the decoherence rate of the the DTLS in free decay given by 

K2 = --Ri2,i2- (55) 

These results accord with those obtained by others.— In general, due to dephasing T2 < 2Ti.— 
The dephasing time can be calculated from Ti and T2 by^ 

i- = i--J-. (56) 

B. Relaxation and decoherence times of the DTLS in the presence of a resonant ac driving 
field 

For an DTLS resonantly driven by an ac field, multiple relaxation and decoherence times are 
required to completely describe time evolution of population and coherence of the qubit.^^ Particu- 
larly, intrinsic and field-induced decoherence times are necessitated to characterize the decoherence. 
In the underdamped regime, the most important regime for the driven qubit, the relaxation time 
Ti is given by 

fi = T-\ (57) 

where, the tilt "~" is used to denote the characteristic times of the resonantly driven qubit and F 
is the relaxation rate of the driven DTLS given by 

r = ^1±^H. (58) 

The intrinsic decoherence time T21 and field-induced decoherence time T22 are given by 

T2i = i^2^ = T2, (59) 

and 

T22 = r-^ = Ti, (60) 

respectively. It is shown that due to the effect of driving fields the characteristic times of the driven 
DTLS are different from those of the DTLS in free decay. The characteristic times of the driven 
DTLS are independent of the field strength in the case of weak driving fields. If the initial state 
of the driven DTLS is an eigenstate such as the ground state, the intrinsic decoherence vanishes 
and the decoherence is completely characterized by the field-induced decoherence. In this case, 
the decoherence rate of the driven DTLS is equal to its relaxation rate as shown by Eq. ()60p . 
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Therefore, the escape rate of quantum phase information from the system equals to the rate of 
energy flowing from the system to the environment. This result is identical with that obtained 
from the non-Markovian approach with inclusion of both bath and qubit dynamicsi^ 
From Eqs. §2\j, ([53]), ([59]), and (l60|) we obtain 



1 _ 1 _ 1 1 _ 3 1 

It shown that T22 = Ti, min(ri,r2) < Ti (t22^ < max(ri,r2), and Ti (^^22^ < min (4Ti/3, 2T<^). 
These results are similar to those obtained by other o^^'^^ and also agrees, within the experimental 
uncertainties, with recent experimental results.— 

VI. DECOHERENCE OF A 2D SQUID FLUX QUBIT COUPLED TO CONTROL AND 
READOUT CIRCUITS 

A number of effects can destroy coherence of a SQUID qubit, of which the fluctuation of 
external circuits used for control and readout of the qubit is one of the important sources of 
decoherence.— In this section, we investigate decoherence of a 2D SQUID flux qubit due to coupling 
with control and readout circuits and effect of driving fields on the decoherence. 

A. Hamiltonian of the microwave-driven 2D SQUID flux qubit 

As shown in FIG. 1, a 2D SQUID flux qubit is a variable barrier rf SQUID in which the single 
Josephson junction in an ordinary rf SQUID is replaced by a low inductance dc SQUIDi^ Suppose 
that the inductance of the superconducting loop of the rf SQUID is L, the critical current is Ic, the 
total magnetic flux enclosed in the rf SQUID loop is ^, the shunt capacitance of each Josephson 
junction in the dc SQUID is C, the critical currents are Id and Ic2, and the total magnetic flux 
enclosed in the dc SQUID loop is ^dc- The Hamiltonian of the 2D SQUID qubit can be written as 

where, nix = 2C$q and ruy = are the masses of the first and second modes, <I>o = h/2e 

and e are the flux quantum and elementary charge, x = ^/^o and y = ^dc/^o are the canonical 
coordinates of the 2D SQUID qubit, px = —ihd/dx and py = —ihd/dy are the canonical momenta 
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conjugate to x and y, and V{x, y) is the potential energy given 



V{x,y) = -f -[x-x,f + '-{y-y,f 



— — cos (27rx) cos (vry) 
47r^ 

+ -^^ sin (27rx) sin (vry) 




(63) 



Here, g = L/21 is the ratio of the inductances of the rf SQUID and dc SQUID, = 2'irLIc/^o, 
5f3i = 2'kL{Ic2 — Id) /^O) and y^ are the fluxes applied to the rf SQUID and dc SQUID in 
the unit of <I>o- The Josephson coupling energy of the rf SQUID is Ej = filc/2e = mu}\(j[3i/A'K'^ , 
where ujlc = l/V LC is the characteristic frequency of the 2D SQUID qubit. The contour of the 
potential energy of the 2D SQUID qubit used in recent experiments^ is plotted in FIG. 2, where 
the parameters of the 2D SQUID qubit are L = 205 pH, C = 32.5 fF, c/ = 17.0, [3^ = 3.7, 5(5 ^ = 0, 
Xe = 0.4991, and ye = 0.387. 

The spectroscopic properties of the 2D SQUID qubit can be obtained by solving the eigenvalue 
equation of Hamiltonian Hs-^^"-^^ In FIG. 3 and FIG. 4, we plot the energy levels and transition 
matrix elements versus the flux applied to the rf SQUID Xg, respectively. When x^, = 0.4991, 
at the position of the arrow, AE'ai = E3 — Ei = 0.259ujlc = 15.95 GHz and IX21/X32I = 0.262, 
where ujlc = 3.874 x 10^^ rad/s. These results are in very good agreement with the experimental 
results.— Note that Xe and ye determine the energy bias and tunnel splitting, respectively. Thus 
the spectroscopic properties of the 2D SQUID qubit can be varied in situ by adjusting Xe and ye- 

In order to manipulate the qubit's states for gate operations, a microwave pulse is applied to 
the 2D SQUID qubit through the first mode. If the interaction of the microwave field and external 
circuits is neglected, the Hamiltonian of the microwave-driven 2D SQUID qubit coupled to the 
external circuits is given by 



where. Hp {x,t) is the interaction of the SQUID qubit and the microwave field and Hj {x,y,t) is 
the interaction of the SQUID qubit and the external circuits (thermal bath). 

If (p (t) is the normalized fiux from the microwave field coupled to the SQUID qubit, which is 
taken to be 



where, cj)^ and are the field strength and frequency, respectively, then Hp {x, t) is given bj*^ 



H (x, t) = Hs {x, y) + Hp (x, t) + Hj {x, y, t) , 



(64) 



Ht) = 4>f, cos {uJf,t) 



(65) 



Hp{x,t) = -^cP[cl)-2{x-Xe)]. 



(66) 
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In general, the external circuits are coupled to the 2D SQUID qubit through both modes. 
However, the coupling through the second mode is negligibly small compared to that through the 
first mode and will be neglected. In the case of weak coupling, the interaction of the 2D SQUID 
qubit and the external circuits is a linear function of x and can be expressed hy^^ 

Hi{x,t) = xq{t), (67) 

where, q (t) is the fluctuation force of the external circuits. If the coupling between the control and 
readout circuits is neglected, q (t) can be decomposed into two parts as 

=?x(t) + ?m(t), (68) 
where, <^x (t) and ?m(i) are the fluctuation forces of the control and readout circuits, respectively. 



B. Spectral densities of the control and readout circuits 

As has been demonstrated, in the case of weak damping, the effect of a thermal bath on 
a quantum system can be characterized by a spectral density^^i^ The spectral density can be 
computed from the frequency-dependent damping coefficient of the quantum Langevin equation 
(see Appendix B). For the 2D SQUID flux qubit, the quantum Langevin equation is given by^ 



where, the damping coefficient Y{t) is the equivalent admittance of the external circuits, V (<^) is 
the potential energy applied to the SQUID qubit, and (t) is the fluctuation force of the external 
circuits. 

By comparison of Eq. ([69]) with Eq. (|Bip . one has q^^,M^C,^^Y/C, and ^ ^ 
Applying these relations to Eq. ()B3p . we obtain the spectral density J {to) of the fluctuation force 

U it) 



J (lo) = hjjYuiuj) 



1 + coth 



(70) 



where, T is the temperature of the external circuits and Yii{uj) is the real part of the frequency- 
dependent equivalent admittance Y(uj) of the external circuits. Y(ij) can be calculated from Y(t) 
by Fourier transform. It can also be calculated directly from the external circuits using circuit 
equations in the frequency domain. For the control and readout circuits given in FIG. 1, Y{uj) has 
been derived (see Appendix C for details) and Yr{lli) is given by 

YR{uj)=YxR{u;) + Y^R{io), (71) 
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where, Yxr{lo) and YmR{uj) are the real parts of the frequency-dependent admittances of the control 
circuit and readout circuit, respectively. They are given by 

and 

(UJ) [1 + (tj) UJ^l 

respectively, where, (to), Gx i^^), F„i {uj), and Gm {'^) are given by Eqs. (jCSp . ()C9p . ()C22p . and 
(IC23l) in Appendix C. Substituting Eq. ^ into Eq. ([TOD, one has 

J (l<j) = Jx (w) + Jm (w) , (74) 

where, Jx {oj) and Jm (w) are the spectral densities of the control circuit and readout circuit, 
respectively. They can be calculated using Eq. ([70]) with Ypi{uj) replaced by Yxr{uj) and YmRii^), 
respectively. 

In FIG. 5 we plot the spectral densities Jx (w), Jm {'^)^ and J {oj) versus the frequency uj for 
the control and readout circuits of the 2D SQUID qubit used in the experiment.^ The parameters 
of the control circuit are Lx = 100 pH, Gx = 25 pF, Rx = 70 Q., and Rxq = 1.0 x 10^ il. The 
parameters of the readout circuit are Lio = L20 = 20 pH, Lji = 100 pH, L j2 = 550 pH, Gm = 20 
pF, Rm = 70 and Rmo = 2.0 x 10^ VL. The mutual inductances Mx = 1.0 pH and Mm = 3.3 
pH, and the temperature T = 30 mK. It is shown that the spectral density J {u)) (also Jx (1^) and 
Jm (i^)) reaches the maximum at lij ~ 1.69liJlc'. At low frequency, the spectral density J {uj) (also 
Jx {u}) and Jra (f^)) approaches a constant value, in particular J (0) 7^ 0. As will be demonstrated, 
the dephasing time T^p is finite if J (0) 7^ 0. Thus the control and readout circuits of the 2D SQUID 
qubit induce both relaxation and dephasing. 

C. Damping rate matrix 

From Eq. (f67|) . the fluctuation force (t) is given by 

dHi{x,t) 1 
= d^ = ~^o^*^- ^^^^ 

This equation provides a relation between the fluctuation force ? (t) and (t) of the external 
circuits. Compared Eq. ([75]) with Eq. ([l3j) . one has A = — l/<I>o. Using Eq. (|l5]) . the spectral 
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density J<j (uj) of the fluctuation force ? (t) can be calculated from the spectral density J (uj) of the 
fluctuation force (t) by 

J,M=$2j(^). (76) 

For the 2D SQUID qubit, Xi = x and X2 = y. Substituting Eqs. dMD and ([TeD into Eq. dH]) we 
obtain the steady damping rate matrix element Rmn m'n' of the 2D SQUID qubit 



*2 

p — 



2/l2 



k 



(77) 



It is a function of transition matrix elements and spectral densities. 

If only considering the two computational levels, the 2D SQUID qubit is equivalent to a DTLS 
and its characteristic (relaxation, decoherence, and dephasing) times can be calculated analytically. 
Substituting Eq. ([TTD into Eqs. ([Ml) and ([55]), then Eqs. ^ and ^ into Eqs. ([52]), dSS]), and 
([56]) . one has 

2 

^r' = ^M'^iJM + JM], (78) 

= ^ + ^ (^11 - ^22)' ^ (0) , (79) 

and 

T^-i = ^(xn-X22)'j(0). (80) 

It is shown that the relaxation and dephasing rates are determined by the spectral densities at 
transition frequency to = |w2i| and low frequency w = 0, respectively . ^^'^'^'^^i^^'^'^i^^ The relaxation 
rate is proportional to the modulus square of transition matrix element |xi2| , while the dephasing 
rate is proportional to the squared difference of average coordinates of the two states (xn — X22)^- 
For a qubit having (xn — X22) = 0, the dephasing is completely suppressed. For the 2D SQUID 
qubit considered here, both J (0) and (xn — ^22)^ are not zero. Thus the control and readout 
circuits will induce phase relaxation. 

For the 2D SQUID flux qubit, Yxr{u)), YmR{uj), and thus Yfi{io) are even functions of uj. 
Substituting Eq. ([70]) into Eqs. ([78]) to ([80]) . the analytical expressions of the characteristic times 
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are simplified to 



= '^huJ2i \X12\^ Yr {uJ2l) coth (^-^^^ , (81) 



= + 4^bT {xn - X22f Yr (0) , (82) 



and 

= -J^s^ (^11 - ^22)' (0) . (83) 

Now the relaxation rate is dominated by the external circuits' admittance at transition frequency 
UJ21 while the dephasing rate by the admittance at w = 0. Furthermore, the dephasing rate is 
proportional to the temperature of the thermal bath. Hence at low temperature the dominating 
source of decoherence is relaxation while at high temperature the main source of decoherence is 
dephasing. These results agree with those obtained by others j^i^^i^ 

Using Eq. (|i9]) . the general form of spontaneous decay rate, we obtain the spontaneous decay 
rate of the 2D SQUID flux qubit 

27r 



^mn — -^RglXmnl {En — Em)YR{uJr. 



1 + coth 



E — E 



(84) 



2kBT 

where, Rq = h/Ae^ is the resistance quantum and ujnm = {En — Em) /h > 0. If the equivalent 
impedance is a resistance i?, then YR{u)nm) = 1/R and the spontaneous decay rate is the same as 
that given by others with different methodj^i^ 



D. Numerical method 



In the Hilbert space spanned by the eigenstates of the qubit's Hamiltonian, the density matrix 
is governed by the master equation (jlSp . This equation can be rewritten in matrix form as 

^ = -iLp + Rp, (85) 

where, p = {p^} is the density matrix, L = + ^^/^'j is the matrix of Liouville superopera- 

tors, R = is the damping rate matrix, fi = mn, and ^' = m'n'. In Eq. (j85|) the Lamb shift 

matrix has been neglected. The matrix L is a time-independent, real, and symmetric matrix. Its 
matrix elements are given by Eqs. (jlSp and (jl9p . The matrix R is also a time-independent and 
real but non-symmetric matrix. Its matrix elements are given by Eq. (I77p . 
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To solve Eq. (j85p . we use the split-operator method. Using this method, the propagation of 
the density matrix is calculated by 

p{t + At) = Pl {t) Pr (t) Pl {t) p (t) , (86) 

where, P^^ (t) and P/j {t) are the propagators. They are given by 

Pe (t) = exp (XQAt) , (87) 

where, A = -i/2 and Q = L (t + At/2) for 9 = L, and A = 1 and Q = R (t + At/2) for Q = R. 

Suppose that Q is an x matrix, and its eigenvalue and right eigenvector are qk and 
Bk = [bik,b2k^ ■ ■ ■ ,b]\ik]^, respectively, where T denotes the transpose of the vector. We assume 
that q = {qkSk'k} and B = {bk'k} are the two matrices constructed respectively by the eigenvalue 
Qk and eigenvector Bk, then the matrix Q can be calculated by 

Q = BqB \ (88) 

where, is the inverse matrix of B. Substituting Eq. (jSSp into Eq. ()87p . the propagator Pe (t) 
is then given by 

Pe (t) = B exp (XqAt) B'^, (89) 

where, exp (AqAt) is a diagonal matrix with non-zero diagonal matrix elements given by 
exp (A(7fcAt). If Q is a Hermitian operator B^^ = B^^ and if Q is a real and symmetric opera- 
tor (e.g., L) = B^. 

E. Effect of driving fields on relaxation and decoherence times 

By numerically solving the master equation in terms of the numerical method introduced in 
the preceding section, population and coherence of the microwave-driven 2D SQUID qubit are 
calculated, from which the relaxation and decoherence times are extracted. 

1. Free decay of the 2D SQUID qubit 

To numerically calculate the relaxation and decoherence times of the 2D SQUID qubit in the 
absence of driving fields (free decay) we assume that the initial state of the qubit is a superposition 
state of the two computational states |1) and |2) with pn (0) = (0) = P21 (0) = P22 (0) = 0-5 
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and Pij{0) = for all the other combination of i and j. To take into account the leakage to non- 
computational states we include four levels (A^ = 4) in the calculation. In the case of weak driving 
fields and weak damping, the calculation with four levels {N = 4) is converged. By numerically 
solving the master equation (j85p with the aforementioned initial state we obtain the population 
and coherence of the 2D SQUID qubit in free decay. Since the coherence is usually a complex and 
fast oscillating, we use I/J12P instead of to estimate the decoherence time. In FIG. 6 and FIG. 
7, we show with the solid lines the evolution of population inversion (p22 ~ Pii) &iid the squared 
modulus of coherence \pi2\ , respectively. 

For the DTLS, the population inversion and squared modulus of coherence in free decay undergo 
simple exponential decays with^ 

P22 - Pii = yi + ^ie-*/"\ (90) 

and 

l/0i2l' = 2/2 + ^2e-2*/-2, (91) 

where, the parameters ti and T2 are the relaxation and decoherence times of the DTLS in free 
decay, respectively. 

To calculate the relaxation and decoherence times of the 2D SQUID qubit in free decay, we fit 
the numerical results of (p22 ~ Pii) \pi2\'^ to the above exponential functions. The results of 
least-square fitting are plotted in FIG. 6 and FIG. 7 with the dashed lines, from which we obtain 
Ti = Ti = 3.429 ps and T2 = T2 = 2.243 ps. The calculated relaxation time is in very good 
agreement with the experimental result Ti = 3.45 ps^ demonstrating the validity of our approach 
and calculation. Using Eq. (j56p . we obtain the dephasing time = 3.333 ps. For comparison, we 
also calculate the relaxation and decoherence times using the analytical expressions of the DTLS 
in free decay given by Eqs. (fSTI) to ([83]) . The results Ti = 3.429 ps, T2 = 2.243 ps, and T^p = 3.333 
ps are exactly same as the numerical results. In addition, since the decoherence time is shorter 
than the relaxation time the dephasing is the main source of decoherence. 

2. Rabi oscillation of the resonantly driven 2D SQUID qubit 

To calculate the relaxation and decoherence times of the 2D SQUID qubit in the presence of a 
microwave field, we assume that the initial state of the qubit is the ground state with p^^ (0) = 1 
and Pjj(O) = for all the other i and j. By numerically solving the master equation (j85p with this 
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initial condition, the population and coherence and thus the relaxation and decoherence times of the 
driven 2D SQUID qubit are calculated. The relaxation and decoherence times of the driven qubit 
depend on the relative value of the field strength to the damping strengthi^ In the underdamped 
regime for which the field strength is larger than the damping strength, the decoherence can be 
decomposed into intrinsic and field-induced ones. If the initial state of the qubit is the ground 
state, the intrinsic decoherence vanishes and the qubit has a single decoherence time which equals 
to the field-induced decoherence time. 

In FIG. 8 and FIG. 9 we plot with the solid lines the evolution of population difference {pn — P22) 
and squared modulus of coherence \pi2\ , respectively, for the 2D SQUID qubit resonantly driven 
by the microwave field with cj)^ = 1.0 x 10^^ and lo^ = u)2i = 0.127ujlc- As shown in these figures, 
both (pii — P22) and \pi2\^ undergo damped Rabi oscillations. 

In the underdamped regime, the population difference and squared modulus of coherence of the 
resonantly driven DTLS from an eigenstate undergo damped Rabi oscillations as^ 

Pii - P22 = yi + zi sin {nt + if^) e-*/^\ (92) 

and 

\P12f = y2 + 22 sin + (/?2) e"*/^2 

+Z3 sin^ [Qt + ip2) e-^*/^2 ^ (-93) 

where, is the Rabi frequency, and ri and T2 are the relaxation and (field-induced) decoherence 
times of the driven DTLS, respectively. 

To extract the relaxation and decoherence times of the driven 2D SQUID qubit, we fit the calcu- 
lated {pii — P22) and 1^12!'^ ^° aforementioned exponentially damped Rabi oscillating functions. 
The results of the best fit are shown in FIG. 8 and FIG. 9 with dashed lines, from which we obtain 
n = 4.016 X IQ-^ujLc, Ti = ?i = 2.689 ^s, and 722 = T2 = 2.682 ^s. 

Using the same procedure, we have calculated the relaxation and decoherence times of the 2D 
SQUID qubit resonantly driven by the microwave fields with different field strengths. The results 
are given m the columns with TV = 4 in TABLE HI To examine the effect of leakage on the relaxation 
and decoherence times, we have also calculated the relaxation and decoherence times only using 
the two computational states. The results are given in the columns with = 2 in TABLE [H In 
the case of weak driving fields and weak damping, the 2D SQUID qubit may be well approximated 
by a DTLS. Form the analytical expressions of relaxation and decoherence times of the driven 
DTLS given by Eqs. and §2^, we obtain fi = T22 = 2.712 ps. They are independent 

of driving field strength. 
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TABLE H] shows that when (j)^ < 1 x the calculated relaxation and decoherence times 

of the driven 2D SQUID qubit are essentially identical. They are independent of the driving 
field strengths, shorter than the relaxation time and longer than the decoherence time of the 2D 
SQUID qubit in free decay. These results accord with those obtained from the calculation with the 
analytical expressions of characteristic times of the driven DTLS. The relaxation and decoherence 
times obtained from the calculation with = 4 are the same as those with N = 2 and those 
with the analytical expressions, demonstrating that both the strong field effect and leakage are 
negligibly small in this case. When i?!*^ > 5 x 10^^, on one hand, the relaxation time obtained from 
the calculation with = 4 equals to that with N = 2, indicating that the leakage does not influence 
the relaxation time. On the other hand, the relaxation time obtained from the calculation with 
A^ = 2 is less than that obtained from the calculation with the analytical expressions, illustrating 
that the strong field effect makes the relaxation time smaller. In contrast, due to the strong field 
effect the decoherence time obtained from the calculation with N = 2 increases slowly with the 
field strength, while due to the leakage the decoherence time obtained from the calculation with 
A^ = 4 decreases with the field strength quickly. In particular, when 0^ > 1 x 10~^ the squared 
modulus of coherence obtained from the calculation with A^ = 4 no longer undergoes the simple 
damped Rabi oscillation as that given by Eq. (j93p . Thus the relaxation time is sensitive to the 
strong field effect while the decoherence time is sensitive to the leakage. 

TABLE I: Numerical results of relaxation and decoherence times (fis) of the 2D SQUID qubit in free and 
driven decays. 



Field strength 


Relaxation 


time 


Decoherence 


time 




A = 4 


A = 2 


A = 4 


A = 2 





3.429 


3.429 


2.243 


2.243 


1 X 10-'^ 


2.712 


2.712 


2.712 


2.712 


5 X 10-^ 


2.712 


2.712 


2.712 


2.712 


1 X IQ-s 


2.712 


2.712 


2.712 


2.712 


5 X 10"^ 


2.706 


2.706 


2.705 


2.713 


1 X 10-^ 


2.689 


2.689 


2.682 


2.716 


5 X IQ-^ 


2.224 


2.224 


1.945 


2.742 


1 X 10-" 


1.480 


1.480 




2.837 
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F. Optimization of the control and readout circuits 

The decoherence of the 2D SQUID qubit strongly depends on the control and readout circuits. 
To optimize the control and readout circuits for long decoherence time, we investigate how the 
characteristic times change with the parameters and temperature of the circuits of the 2D SQUID 
qubit in free and driven decays. The strength of the resonant microwave field used is cj)^ = Ix 10~^. 

1. Characteristic times versus mutual inductances 

In FIG. 10 we show the characteristic times versus the mutual inductance between the 2D 
SQUID qubit and the control circuit. In this figure, Ti and T2 are the relaxation and decoherence 
times of the qubit in free decay, and Ti and T22 are those of the qubit in driven decay. It is shown 
that for all the values of Mx shown in this figure Ti ~ T22 and min(Ti,T2) ^ Ti < max(ri,T2). 
These results agree with the predictions from the analytical expressions for the DTLS given by Eq. 
(j6ip and hold, as will be shown, for different control- and readout-circuit parameters. It is also 
shown that when is less than 0.6 pH T2 > Ti, when is larger than 0.6 pH Ti > T2, and 
when Mx is equal to 0.6 pH Ti ~ T2. Hence the relaxation is the dominating source of decoherence 
for smaller Mx while the dephasing is the main source of decoherence for larger Mx- In addition, 
when Mx < 0.1 pH the characteristic times do not change with Mx and when Mx > 0.1 pH the 
characteristic times decrease with Mx monotonically. These results can be analyzed by using the 
analytical expressions of characteristic times for the DTLS. For example, from Eqs. (jSip . (j72p . 
(jCSp . and ()C9p . for small and large Mx, the relaxation time Ti of the qubit in free decay can be 
well approximated by 

Tf 1 ~ ao (1 + boM^) , (94) 

where qq and 60 are two parameters independent of Mx- Eq. ()94p shows that Ti is a constant for 
Mx < l/Vbo and Ti oc M"^ for Mx > l/Vbo, which are in good agreement with the numerical 
results shown in FIG. 10. The changes of T2, Ti, and T22 with Mx are very similar to that of Ti 
with Mx- 

In FIG. 11 the characteristic times are shown versus the mutual inductance M^ between the 
2D SQUID qubit and the readout circuit. It is shown that when Mm is less than 6 pH Ti > T2 
and the dephasing is the main source of decoherence, when M^ is larger than 6 pH T2 > Ti and 
the relaxation is the dominating source of decoherence, and when ~ 6 pH Ti — T2. It is 
also shown that when M^, < 0.5 pH the characteristic times do not change with M^, and when 
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Mm > 0.5 pH the characteristic times decrease monotonically with Mm- These behaviors are very 
similar to those of the characteristic times versus M^. Thus the dependence of the characteristic 
times on Mm can also be well approximated by an equation analogous to Eq. (p^ . 

2. Characteristic times versus readout- circuit parameters 

In FIG. 12 to FIG. 15, we plot the characteristic times versus the inductances Lji, Lj2, Lio, 
and L20 of the readout circuit, respectively. It is shown that the characteristic times change 
dramatically with Lji, Lj2, and Lio but decrease slowly with L2o- When Lji {Lj2 or Lio) is less 
than 550 pH (100 pH or 470 pH) the characteristic times increase with Lji {Lj2 or Liq). When 
Lji {Lj2 or Liq) equals to 550 pH (100 pH or 470 pH) the characteristic times reach the maxima. 
After Lji {Lj2 or Liq) is larger than 550 pH (100 pH or 470 pH) the characteristic times decrease 
with Lji {Lj2 or Lio). In contrast, the characteristic times do not change much with L20 for very 
small and very large L20 and decrease with L20 for moderately large L2Q. In addition, for all the 
values of inductances shown in these figures, Ti > T2. Thus the dephasing is the dominating source 
of decoherence 

To gain insights into mechanisms behind these behaviors, we analyze these results using the 
analytical expressions of characteristic times for the DTLS. From Eqs. (jSip . (j73p . and ()C22p . at 
the adjacency of AL = the relaxation time Ti of the qubit in free decay can be approximated by 

T^-^ ~ ai + 61 (AL)^ , (95) 

where, AL = (L20 + Lj2) — {Liq + Lji), ai is a parameter independent of AL and bi is a parameter 
slowly and smoothly varying with AL. Eq. (j95p shows that Ti decreases with AL and reaches the 
maximum at AL = 0. 

For the results in FIG. 12, the position of the maxima of the characteristic times is at Lji = 550 
pH. This result is in very good agreement with the prediction from Eq. (|95p since when Lji = 550 
pH AL = 0. For the results in FIG. 13 and FIG. 14, we also have AL = at the positions of the 
maxima of the characteristic times. As for the results in FIG. 15, AL is always larger than 430 pH 
since AL = L20 + 430 pH and L20 > 0. Thus the characteristic times do not have maxima on L2q. 

In FIG. 16 we exhibit the characteristic times versus the capacitance Cm of the readout circuit. 
When Cm < 0.01 pF the characteristic times do not change with Cm, when 0.01 pF < Cm < 2 pF 
the characteristic times decrease with Cm , and when Cm > 2 pF the characteristic times tend to 
constants. For the values of Cm shown in this figure, T2 < Ti and thus the dephasing is the main 
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source of decoherence. 

In FIG. 17 and FIG. 18, we show the characteristic times versus the resistances Rm and Rmo 
of the readout circuit, respectively. FIG. 17 shows that the characteristic times increase with Rm 
when Rm < 1000 Q and do not change with Rm when Rm > 1000 J7. When Rm is less than 21 

T2 > Ti and the relaxation is the main source of decoherence, when Rm is larger than 21 Q 
Ti > T2 and the dephasing is the dominating source of decoherence, and when Rm equals to 21 Q 
Ti = T2. FIG. 18 shows that the characteristic times increase with Rmo when Rmo < 2000 Q, and 
tend to constants after Rmo > 2000 U. For the values of RmO shown in the figure, Ti > T2 and the 
dephasing is the main source of decoherence. 

3. Characteristic times versus control- circuit parameters 

In FIG. 19 we plot the characteristic times versus the inductance of the control circuit. It is 
shown that when < Ix 10^ pH the characteristic times are constants, when 1 x 10^ < L^, < 1 x 10^ 
pH the characteristic times increase with L-x, and when Lj; > 1 x 10^ pH the characteristic times 
do not change with L^. For all the values of shown in this figure T2 < Ti, demonstrating that 
the dephasing is the main source of decoherence. 

In FIG. 20 the characteristic times are plotted versus the capacitance Cx of the control circuit. 
It is shown that the changes of the characteristic times with Cx are very similar to those with Cm- 
When Cx < 0.05 pF the characteristic times do not change with Cx, when 0.05 < Cx < 2.5 pF 
the characteristic times decrease with Cx, and when Cx > 2.5 pF the characteristic times tend to 
constants. For all the values of Cx shown in this figure T2 < Ti and thus the dephasing is the main 
source of decoherence. 

In FIG. 21 and FIG. 22 the characteristic times are plotted versus the resistances Rx and Rxq of 
the control circuit, respectively. It is shown from FIG. 21 that when Rx is less than 20 Q T2 > Ti 
and the relaxation is the main source of decoherence, when Rx is larger than 20 f2 Ti > T2 and 
the dephasing is the dominating source of decoherence, and when Rx equals to 20 Ti ~ T2. For 
Rx < 2000 O the characteristic times increase with Rx and for Rx > 2000 O the characteristic 
times tend to constants. The changes of the characteristic times with Rxo shown in FIG. 22 are 
very similar to those with Rx- When Rxo is less than 2200 $7 Ti > T2 and the dephasing is the 
main source of decoherence, when Rxo is larger than 2200 0, T2 > Ti and the relaxation is the 
dominating source of decoherence, and when Rxo is equal to 2200 Q Ti = T2. For Rxo < 3 x 10^ Q, 
the characteristic times increase with Rxo while for Rxo > 3 x lO'' Q, the characteristic times tend 



29 



to constants. 



4- Characteristic times versus temperature 

Finally, we plot the characteristic times versus the temperature T of the external circuits in 
FIG. 23. It is shown that in general the characteristic times decrease with T. When T ^ the 
characteristic times tend to constants and when T ^ oo the characteristic times are inversely 
proportional to T which is the results of classical mechanics. These results agree with the predic- 
tions from the analytical expressions of characteristic times for the DTLS. It is also shown that 
at the lower temperature when T < 10.3 mK T2 > Ti and the relaxation is the main source of 
decoherence, while at the higher temperature when T > 10.3 mK Ti > T2 and the dephasing is the 
dominating source of decoherence. These results also agree with the predictions from the analytical 
expressions for the DTLS^ and with those obtained by others.— i^^i^. 



VII. CONCLUSION 



In summary, to investigate the environment-induced decoherence in realistic gate operations 
of solid-state qubits, we present a general theory for the treatment of decoherence of a multilevel 
quantum system of many degrees of freedom interacting with a multibath reservoir and driven by 
ac fields. In this theory, the system is described by the reduced density operator governed by the 
master equation. The effect of the environment on the system is characterized by the spectral 
density through the dissipation superoperator. The effects of driving field and leakage due to the 
coupling with both the driving field and environment are included in this theory. In the Hilbert 
space spanned by the eigenstates of the system's Hamiltonian, the reduced density operator is 
represented by the density matrix and the dissipation superoperator by the dissipation matrix. 
The diagonal and off-diagonal matrix elements of the density matrix stand for the population 
and coherence of the system, respectively. The dissipation matrix can be decomposed into Lamb 
shift matrix and damping rate matrix. They are determined by the transition matrix elements 
of the system and the spectral density of the environment. The effect of the Lamb shift matrix 
on the system is analogous to an extra field. In the case of weak damping the Lamb shift matrix 
is extremely small compared to the driving field and are neglected. In the study of decoherence 
of a qubit, for which the long-time behavior of the qubit is significant, the damping rate matrix 
is replaced by a steady one. For the thermal bath, the spontaneous decay rate and stimulated 
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transition rate are derived from the damping rate matrix. They accord with those obtained by the 
others with different methods and obey the detailed balance principle. 

For an DTLS, the characteristic times in free and resonantly driven decays are expressed by 
the analytical expressions. The decoherence of the driven qubit can be decomposed into intrinsic 
and field- induced ones. The intrinsic decoherence time equals to the decoherence time of the qubit 
in free decay and the field-induced decoherence time equals to the relaxation time of the driven 
qubit. In the case of weak driving fields, the relaxation and thus the field-induced decoherence 
times of the driven qubit are independent of the field strengths and they are always in between the 
relaxation and decoherence times of the qubit in free decay. 

For demonstration, we have applied the dissipative theory to simulate the dissipation process 
of the 2D SQUID qubit coupled to the external circuits in free decay. The energy levels, transition 
matrix elements, and relaxation time are in very good agreement with the experimental results. 
Wc have also applied the dissipative theory to investigate the effect of driving field and leakage on 
the decoherence of the 2D SQUID qubit coupled to the external circuits and resonantly driven by 
the microwave field. In the case of weak driving fields, the relaxation and decoherence times of the 
driven qubit are identical. They are independent of the driving field strength and in between the 
relaxation and decoherence times of the qubit in free decay. These results agree with the analytical 
results obtained from the analytical expressions of characteristic times for the DTLS. In the case of 
a little bit stronger driving fields, the relaxation time is sensitive to the strong field effect while the 
decoherence time is sensitive to the leakage. In addition, for the qubit in free decay, the relaxation 
is the main source of decoherence at the low temperature while the dephasing is the dominating 
source of decoherence at the high temperature. 

To optimize the external circuits for long decoherence time, we have investigated the charac- 
teristic times of the 2D SQUID qubit change with the parameters and temperature of the control 
and readout circuits. We found that the characteristic times decrease with the mutual inductances 
and capacitances, increase with the resistances, and change dramatically with the inductances, in 
particular, the inductances of the readout circuit. To gain longer decoherence time, the coupling 
of the 2D SQUID qubit and the external circuits should be weak, the capacitances of the external 
circuits should be smaller, the resistances should be larger, the inductance of the control circuit 
should be properly larger, and in particular, to reduce the damping due to the readout circuit the 
total inductance of the left branch of the readout circuit should balance with that of the right 
branch. 
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APPENDIX A: MASTER EQUATION OF AN OPEN QUANTUM SYSTEM 



To derive the master equation for the reduced density operator of an open quantum system, we 
work in interaction picture. In this picture the density operator of the global system rj{t) and the 
reduced density operator of the quantum system are defined by 



7y(t) =exp [i iCs + CR)t] 7]{t), 



(Al) 



and 



p{t) = exp [iCst]p{t), 



(A2) 



where, rj{t) and p{t) are the density operator of the global system and the reduced density operator 
of the quantum system in the Schrodinger picture, respectively. The relation between and r]{t) 
is also given by Eq. Applying Eq. ()Aip to Eq. ([2]) we obtain the Liouville-von Neumann for 
rj(t) in the interaction picture 



dt 



CF{t)+Ci{t) ?,{t), 



(A3) 



where, Cp (t) and £/ (t) are the presentations of Cp and jOj in the interaction picture. They are 
defined by 



Cpit) = exp {iCst) CpQ'x.V {—i^st) > 



(A4) 



and 



£/(t) = exp [i (£5 + Cr) t] Ci exp [-i (Cs + Cr) t] 



(A5) 



respectively. 

The formal solution of Eq. (jA3|) is 



r](t) = T exp 



-i (CFir) + CiiT))dr 



(A6) 
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where, T is the time-ordering operator and ^^(O) is the initial density operator of the global system 
which is equal to ??(0) from Eq. ()A1|) . 

Suppose the reservoir is uncorrelated with the system at i = 0. In this case, the initial density 
operator of the global system r/ (0) can be written as^ 

r/(0) = p(0)(T(i?), (A7) 

where, p (0) is the reduced density operator of the system at t = and a (R) is the density 
operator of the reservoir. For thermal baths at temperature T, for example, a (R) is determined 
by the Boltzmann distribution^ 

where, /c_b is the Boltzmann constant. We also assume that (1) both the interactions Hj and Hp 
are weak so that the Born perturbation approximation can be applied to the exponential function 
of Eq. (|A6p :— (2) the reservoir is sufficiently large and its states are unperturbed by the coupling 
with the system and obey Gaussian statistics so that a (R) is time-independent, Tir [a (R)] = 1, 
and T'rii[Hj{t)(T (R)] = 0; and (3) the characteristic time of correlation is much less than the 
relaxation time of system so that the change of p{t) is slow and p(r) can be replaced by p{t) in 
the integral over the correlation time.— Under the these assumptions the system will perform a 
Markovian process. Substituting Eq. (jA6j) into Eq. (jll), calculating the trace over the reservoir 
by means of the cumulant expansion method^^ under the aforementioned assumptions, and using 
the initial condition of Eq. ()A7p . we finally obtain an equation for in a series of cumulant, 
which is given, up to the second order of cumulant, by 

=f exp|-i^ CF{T) + iVi{T) dT^p{0), (A9) 

where, 'Dj{t) is the dissipation superoperator in the interaction picture given by 

Vjit) = - f Ttr \Ci (t) Ci (r) a (R)] dr. (AlO) 

Jo L J 

If making transformations rj{t) 'p{t) and Ci{t) — > iDi{t) in Eq. ()A6P one can get Eq. ()A9p . 
Thus if making the same transformations in Eq. ()A3p one obtains an equation of motion for 'p{t) 

^ = -iCpimt) + vi{t)~p{t\ (All) 

This equation is recognized the master equation in the interaction picture. 
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Applying Eq. ()A2p to Eq. (jAlip . one obtains the master equation for p{t) in the Schrodinger 
picture 

^ = -z [Cs + Cpit)] p{t) + Vj{t)p{t), (A12) 

where, Vi (t) is the representation of T)i[t) in the Schrodinger picture given by 

Vi{t) = - f Ttr [CI (t, r) £f (r, t) a {R)] dr, (A13) 
Jo 

with 

£f (ii, ts) = exp ^/ ih) exp [^£5*2] • (A14) 

APPENDIX B: QUANTUM LANGEVIN EQUATION AND FLUCTUATION- 
DISSIPATION THEOREM 

We consider a quantum particle of mass M moving in a potential V {q) and linearly coupled to 
a thermal bath at temperature T. In phase space, the motion of the quantum particle is described 
by the quantum Langevin equatioii^ 

-^-^/>-<'-^)^-^-^«> 

where, q{t) is the coordinate operator of the particle, 7 (t) is the damping coefficient, and ^ (t) 
is the fluctuation force of the thermal bath. The quantum Langevin equation provides a relation 
between the damping coefficient 7 (t) and the fluctuation force ^ (t) of the thermal bath. 

In frequency domain, the frequency-dependent damping coefficient 7 (lo) is calculated from 7 (t) 
by Fourier transform 

/+00 
7 (t) exp {-iojt) dt. (B2) 
-00 

The spectral density of the fluctuation force ^ (t) at temperature T, (uj), can be calculated from 
7^ (uj) in terms of the quantum fluctuation-dissipation theorem^i^ 



(1^) = MHuj^j^ {uj) 



1 + coth 



(B3) 



At high temperature limit, T ^ fuv/ks and (w) in Eq. ()B3P is simplified to the result of the 
quasiclassical and classical fluctuation-dissipation theorem^ 

(uj) = Mhwjji {uj) coth 

~ 2MkBT-iji{uj) . (B4) 
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It is shown that for a classical particle {to) is proportional to the temperature. In contrast, at 
low temperature limit, T <C huj/ks and (a;) in Eq. (|B3p is simplified to 

[uj) = Mfuvyji (w) , (B5) 

which is independent of the temperature. 

The autocorrelation function of the fluctuation force ^ (t), (t) = (t) ^ (0)), is calculated 
from (uj) by using Eq. (|12p . If jf, (uj) is an even function of uj, the autocorrelation function of 
a quantum particle is given hy^ 



vr 



hoo 



hM . 

Ji {t) = / ^^iR M [i sm [uJt) 







+ coth ( 7-^^ ) cos (ujt)] dio. (B6) 
\2kBT J 

And the autocorrelation function of a quasiclassical particle is given by 

X COS (iot) duj. (B7) 

From Eq. (|B7p the autocorrelation function of the quasiclassical particle is real and reversible, i.e., 
(0^ (0)) = i—t)C (0)). However, from Eq. ()B6p the autocorrelation function of the quantum 
particle is complex and irreversible. 

APPENDIX C: SPECTRAL DENSITY OF THE CONTROL AND READOUT CIR- 
CUITS 

For the 2D SQUID flux qubit shown in FIG. 1, the thermal bath consists of the control and 
readout circuits. They produce thermal noises for the 2D SQUID qubit. 

1. Spectral density of the control circuit 

As shown in FIG. 24, the control circuit is coupled to the rf SQUID through the loop L with 
mutual inductance Mx- The equivalent impedance Zi of the parallel branches CxRx and Rxo is 
given by 

z, = Rmu;)+ . I (CI) 

jujCx^ (uj) 



35 



where, (lu) is the frequency-dependent equivalent resistance given by 

" ^ ^ l + a;2C|(i2, + i?,of ^ ^ 

and Cx'^ (lj) is the frequency-dependent equivalent capacitance given by 

When Rxo ^ oo, iJl'^ = Rx and Cj'^ = Cx- When a; ^ 0, i?^'! = Rxo and w^^q = l/Ca^i^^Q. 
The circuit equations of the control circuit are 

jcoLI - jujMxIo = U, (C4) 

- juMxI + ijuLx + Zi)Io = 0. (C5) 
Prom these equations, the equivalent impedance Zx {co) of the control circuit is given by 

TT . .2 /1^2 

Zx {^) = y= ju^L + . , % ■ (C6) 

Prom this equation, the frequency-dependent equivalent admittance Yx{u) = 1/Zx{lo) is obtained 
and its real part Yxr (co) is given by 

YM^) = ./r!^\.. , (C7) 



where. 



and 



When Rxo — oo 



io^ + Gx{ioy 



F.H = ^^, (C8) 
or r2 / 1 \ 



Vx = Ml - LLx. (CIO) 



M^.R.. 



Fxiu;) = -^, (Cll) 

vi 



and 



2L L2 1 



When oj ^ 



and 



F.(-) = ^^, (C13) 
^5 



Gxiu) = ^^. (C14) 
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2. Spectral density of the readout circuit 



As shown in FIG. 25, the readout circuit is coupled to the rf SQUID through the loop L with 
mutual inductance M^. The equivalent impedance Z-2 of the parallel branches CmRm and RmG is 
given by 

1 



Z2 = W^{U) + 



juC^ {uj) ' 

where, Rm {^) is the frequency-dependent equivalent resistance given by 

{Rm + -Rmo) r> 

and Cm i^^) is the frequency-dependent equivalent capacitance given by 



(C15) 



(C16) 



1 + UJ {Rm + RmoY 

Lv'^CmRmO 



(C17) 



When Rmo oo, R^ = Rm and = Cm- When w ^ 0, R^ = Rmo and lo'^C^ = l/CmRmo- 
The circuit equations of the readout circuit are given by 



juLI + juj^-h - joj—^h = U, 



(C18) 



and 



■ T T ^ ■ ■ T T ■ Mmj 

juLiIi + J^—^I = - JLO—^I, 



(C19) 



(C20) 



where, Li = Lio -|- Lji and L2 = L20 + Lj2. Prom these equations the equivalent admittance of 
the readout circuit Ym {<jj) = I/U can be obtained. Its real part is given by 

1 



YmR (t^) 



where, Fm {uj) and Cm {i^) are given by 



Fm {co) [1 + Cm (w) io' 



21 ' 



(C21) 



Fm{u;) = R'^ 



Mrr 



AL 



(C22) 



and 



Gm (t^) 



rjeq2 



Ln l-k\ 



i) 



1 - kl ^'C^r^t 



(C23) 



37 



where, AL = L2 - Li, L^c = Li + L2, L\\ = LiLi/Loc, kj^ = M^/LLoc, and kj = M^/ALL\\. 
When Rjjio — 00 



and 



When a; — 



and 




(C24) 



(C25) 



(C26) 



(C27) 
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Figure Captions 

FIG. 1 Sketch of the 2D SQUID qubit inductively coupled to the control and readout circuits. 

FIG. 2 (Online color) Contour of the potential energy of the 2D SQUID qubit. 

FIG. 3 (Online color) Energy levels of the 2D SQUID qubit versus Xe- 

FIG. 4 (Online color) Transition matrix elements of the 2D SQUID qubit versus Xg- 

FIG. 5 (Online color) Jxi^^), Jm{i^), and J (to) versus uj for the external circuits of the 2D SQUID 
flux qubit at T = 30 mK. 

FIG. 6 (Online color) Evolution of population inversion of the 2D SQUID qubit in free decay. 
The solid and dashed lines are the numerical and fitting results, respectively. 

FIG. 7 (Online color) Same as FIG. 6 but for the squared modulus of coherence. 

FIG. 8 (Online color) Evolution of population difference of the 2D SQUID qubit resonantly 
driven by the microwave field with 0^ = 1.0 x 10~^ and a;^ = 0^21 = 0.127a; lc. The solid and 
dashed lines are the numerical and fitting results, respectively. 

FIG. 9 (Online color) Same as FIG. 8 but for the squared modulus of coherence. 

FIG. 10 (Online color) Characteristic times versus the mutual inductance between the 2D 
SQUID qubit and the control circuit. In this figure, Ti and T2 are the relaxation and decoherence 
times of the qubit in free decay, and Ti and T22 are the relaxation and decoherence times of the 
qubit in driven decay. 

FIG. 11 (Online color) Same as FIG. 10 but for characteristic times versus the mutual inductance 
Mjn between the 2D SQUID qubit and the readout circuit. 

FIG. 12 (Online color) Same as FIG. 10 but for characteristic times versus the inductance Lji 
of the first junction of the readout circuit. 

FIG. 13 (Online color) Same as FIG. 10 but for characteristic times versus the inductance Lj2 
of the second junction of the readout circuit. 

FIG. 14 (Online color) Same as FIG. 10 but for characteristic times versus the inductance Lio 
of the readout circuit. 

FIG. 15 (Online color) Same as FIG. 10 but for characteristic times versus the inductance L20 
of the readout circuit. 

FIG. 16 (Online color) Same as FIG. 10 but for characteristic times versus the capacitance Cm 
of the readout circuit. 

FIG. 17 (Online color) Same as FIG. 10 but for characteristic times versus the resistance Rm 
of the readout circuit. 
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FIG. 18 (Online color) Same as FIG. 10 but for characteristic times versus the resistance Rmo 
of the readout circuit. 

FIG. 19 (Online color) Same as FIG. 10 but for characteristic times versus the inductance Lx 
of the control circuit. 

FIG. 20 (Online color) Same as FIG. 10 but for characteristic times versus the capacitance Cx 
of the control circuit. 

FIG. 21 (Online color) Same as FIG. 10 but for characteristic times versus the resistance Rx of 
the control circuit. 

FIG. 22 (Online color) Same as FIG. 10 but for characteristic times versus the resistance Rxo 
of the control circuit. 

FIG. 23 (Online color) Same as FIG. 10 but for characteristic times versus the temperature T. 
FIG. 24 The control circuit of the 2D SQUID flux qubit and its equivalent admittance Yx{uj). 
FIG. 25 The readout circuit of the 2D SQUID flux qubit and its equivalent admittance l^(u;). 
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